# This file uses border_tiles.R to generate figures 8-11, Appendix H figures, Appendix J figures

rm(list=ls())

library(openxlsx)
library(ggplot2)
library(ggthemes)
library(reshape2)
library(tokenizers)
library(dplyr)
library(tidyverse)
library(countrycode)
library(stringr)
library(survival)
library(zoo)
library(readr)
library(readstata13)
library(ggthemes)
library(tidyr)
library(sp)
library(plyr)
library(foreign)
library(bife)
library(estprod)

#load country-year analysis data
load("border_tiles.RData")

###########
# SETUP #
###########
border_tiles$Net_Migration_Prop <- border_tiles$Net_Migration/border_tiles$Population

# logging
border_tiles$Log_Population <- log(border_tiles$Population)
border_tiles$Log_GDP <- log(border_tiles$GDP)

# rescaling
border_tiles$Political_Globalization <- border_tiles$Political_Globalization/100
border_tiles$Economic_Globalization <- border_tiles$Economic_Globalization/100

# add a version of democracy with missing coded to non-democracy - note that this changes the interpretation a bit
border_tiles$Democracy <- ifelse(border_tiles$Polity < 6 | is.na(border_tiles$Polity), 0, 1)

# binarized
border_tiles$mention_binary <- ifelse(border_tiles$tiles_borders > 0, 1, 0)

# binarized, undergrad
border_tiles$mention_binary_undergrad <- ifelse(border_tiles$tiles_borders_undergrads > 0, 1, 0)

# lags
border_tiles <- border_tiles[order(border_tiles$country, border_tiles$year), ]

border_tiles$lagged_mention_binary <- panel_lag(border_tiles$mention_binary, id=border_tiles$country, time=border_tiles$year, lag=1)

border_tiles$lagged_pop <- panel_lag(border_tiles$Population, id=border_tiles$country, time=border_tiles$year, lag=1)
border_tiles$lagged_gdp <- panel_lag(border_tiles$GDP, id=border_tiles$country, time=border_tiles$year, lag=1)
border_tiles$lagged_eg <- panel_lag(border_tiles$Economic_Globalization, id=border_tiles$country, time=border_tiles$year, lag=1)
border_tiles$lagged_pg <- panel_lag(border_tiles$Political_Globalization, id=border_tiles$country, time=border_tiles$year, lag=1)

# net migration
border_tiles <- border_tiles %>% group_by(country) %>% dplyr::mutate(Net_Migration_Imputed = na.approx(Net_Migration, na.rm=FALSE, method='linear'))
border_tiles <- border_tiles %>% group_by(country) %>% dplyr::mutate(Net_Migration_Prop_Imputed = na.approx(Net_Migration_Prop, na.rm=FALSE, method='linear'))

border_tiles$Lagged_Net_Migration <- panel_lag(border_tiles$Net_Migration_Imputed, id=border_tiles$country, time=border_tiles$year, lag=1)


# neighbor data
load('binary_distance.RData')
dyads <- data_binary
rm(data_binary)

# adding disputes
border_tiles <- dyads %>% filter(conttype == 1) %>% left_join(border_tiles %>% select(country, year, Civil_Disputes), 
                                                              by=c('country_2' = 'country', 'year')) %>%
  select(country_1, year, Civil_Disputes) %>% 
  group_by(country_1, year) %>% 
  dplyr::summarise(Civil_Disputes_Neighbor = sum(Civil_Disputes, na.rm=TRUE)) %>% 
  right_join(border_tiles, by=c('country_1'='country', 'year')) %>%
  dplyr::rename(country = country_1)

#####################
# MENTIONS MODELING #
# Figure 8, H1, 9   #
#####################
library(lfe)
library(fixest)

lfe_mod_mentions <- felm(mention_binary_undergrad ~ Interstate_Disputes + 
                           Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                           I(10*Economic_Globalization) + I(10*Political_Globalization) +
                           I(Net_Migration_Prop_Imputed*100)*Polity + Civil_Disputes_Neighbor | country + year | 0 | country + year, data = border_tiles)

summary(lfe_mod_mentions)

s <- summary(lfe_mod_mentions)$coefficients
coeffs_mentions <- data.frame(s)
coeffs_mentions$varname <- rownames(s)
coeffs_mentions$model <- '0-year'

# lagged versions
border_tiles <- border_tiles[order(border_tiles$country, border_tiles$year), ]
for(k in 1:3){
  border_tiles_lagged <- border_tiles %>% select(Interstate_Disputes,
                                                 Civil_Disputes_Neighbor,
                                                 Log_Population,
                                                 Log_GDP,
                                                 Economic_Globalization,
                                                 Political_Globalization,
                                                 Net_Migration_Prop_Imputed,
                                                 Polity,
                                                 mention_binary_undergrad, 
                                                 year, country) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, mention_binary_undergrad), .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_mentions <- felm(mention_binary_undergrad ~ Interstate_Disputes + 
                             Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                             I(10*Economic_Globalization) + I(10*Political_Globalization) +
                             I(Net_Migration_Prop_Imputed*100)*Polity + Civil_Disputes_Neighbor | country + year | 0 | country + year, 
                           data = border_tiles_lagged)
  
  # comparisons between US/NOR migration and interstate disputes
  fixest_mod_mentions <- feols(mention_binary_undergrad ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                                 Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                                 Economic_Globalization + Political_Globalization +
                                 Net_Migration_Prop_Imputed*Polity | country + year, data = border_tiles_lagged)
  
  print(k)
  us_pred <- predict(fixest_mod_mentions, 
                     border_tiles_lagged %>% filter(country == 'USA', year == 2016))
  nor_pred <- predict(fixest_mod_mentions, border_tiles_lagged %>% filter(country == 'NOR', year == 2016))
  
  us_pred_new <- predict(fixest_mod_mentions,
                         border_tiles_lagged %>% filter(country == 'USA') %>%
                           mutate(Net_Migration_Prop_Imputed = Net_Migration_Prop_Imputed + sd(Net_Migration_Prop_Imputed, na.rm=TRUE)) %>%
                           filter(year== 2016))
  
  nor_pred_new <- predict(fixest_mod_mentions,
                          border_tiles_lagged %>% filter(country == 'NOR') %>%
                            mutate(Net_Migration_Prop_Imputed = Net_Migration_Prop_Imputed + sd(Net_Migration_Prop_Imputed, na.rm=TRUE)) %>%
                            filter(year== 2016))
  
  us_pred_inter <- predict(fixest_mod_mentions,
                           border_tiles_lagged %>% filter(country == 'USA') %>%
                             mutate(Interstate_Disputes = Interstate_Disputes + 1) %>%
                             filter(year== 2016))
  
  nor_pred_inter <- predict(fixest_mod_mentions,
                            border_tiles_lagged %>% filter(country == 'NOR') %>%
                              mutate(Interstate_Disputes = Interstate_Disputes + 1) %>%
                              filter(year== 2016))
  
  print((us_pred_new - us_pred)/(us_pred_inter - us_pred))
  print((nor_pred_new - nor_pred)/(nor_pred_inter - nor_pred))
  
  
  s <- data.frame(summary(lfe_mod_mentions)$coefficients)
  s$model <- glue::glue('{k}-year')
  s$varname <- rownames(s)
  coeffs_mentions <- rbind(coeffs_mentions, s)
}

coeffs_mentions <- coeffs_mentions %>% mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, model) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          `I(Net_Migration_Prop_Imputed * 100)` = 'Net Migration (%)',
                          `I(Net_Migration_Prop_Imputed * 100):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(varname = factor(varname, levels=rev(c('Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Economic Globalization', 'Political Globalization', 
                                                'log(Population)', 'log(GDP/Capita)'))))


# fig 8
gg <- ggplot(coeffs_mentions %>% filter(!varname %in% c('log(Population)', 'log(GDP/Capita)'), model == '1-year')) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper), width=0.1, position=position_dodge(width=0.3)) + 
  geom_point(aes(y=varname, x=Estimate, group = model), position = position_dodge(width=0.3), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  # ggtitle('Mention DV, country-year (LPM, country/year FE)') + 
  guides(color = guide_legend(reverse = TRUE)) +
  scale_color_discrete(name = 'Lags')

ggsave(gg, filename = 'Figure 8.png', dpi=1000, width = 5, height=4, bg="white")

# Figure H1: Full results for in-text mention model (Text Figure 8)
gg <- ggplot(coeffs_mentions) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=model), width=0.3, position=position_dodge(width=0.3)) + 
  geom_point(aes(y=varname, x=Estimate, group = model), position = position_dodge(width=0.3), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  # ggtitle('Mention DV, country-year (LPM, country/year FE)') + 
  guides(color = guide_legend(reverse = TRUE)) +
  scale_color_discrete(name = 'Lags')

ggsave(gg, filename = 'Figure H1.png', dpi=1000, width = 5, height=4, bg="white")

# fig 9 interaction plot (net migration, one year lagged model)
lower_ci <- function(x){0.0044496     + x*0.0011278     - 1.96*sqrt(1.120306e-05 + (x^2)*1.555561e-07 + 2*x*1.008158e-06)}
upper_ci <- function(x){0.0044496     + x*0.0011278     + 1.96*sqrt(1.120306e-05 + (x^2)*1.555561e-07 + 2*x*1.008158e-06)}

interaction_df <- data.frame('Polity' = -10:10, 'Effect_Net_Migration' = 0.0044496 + (-10:10)*0.0011278,
                             'Upper_CI' = upper_ci(-10:10), 'Lower_CI' = lower_ci(-10:10))

gg <- ggplot(interaction_df) + geom_point(aes(x=Polity, y=Effect_Net_Migration)) + geom_errorbar(aes(x=Polity, ymin=Lower_CI, ymax=Upper_CI)) + 
  geom_hline(aes(yintercept=0), linetype='dotted') + theme_minimal() + 
  xlab('Polity') + ylab('Estimated Net Migration Effect') + 
  theme(panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
        panel.background = element_blank(), axis.line = element_line(colour = "black")) +
  geom_rug(data=border_tiles, aes(x=Polity, y=0), alpha=0.01, sides='b', position='jitter') +
  ylim(-0.015, 0.03)

ggsave(gg, filename = 'Figure 9.png', dpi=1000, width = 5, height=4, bg="white")

######################
# SENTIMENT MODELING # 
# Figure 10, H2, 11, J
######################
# world sentiment
border_tiles <- border_tiles %>% group_by(year) %>% 
  dplyr::summarise(world_sentiment = mean(mean_predicted_sentiment, na.rm=TRUE)) %>%
  right_join(border_tiles, by='year')

border_tiles <- border_tiles %>% group_by(year) %>% 
  dplyr::summarise(world_anxiety = mean(mean_anxiety, na.rm=TRUE)) %>%
  right_join(border_tiles, by='year')


# region sentiment/anxiety
border_tiles$region <- countrycode(border_tiles$country, 'iso3c', 'un.regionsub.name')
border_tiles$region[border_tiles$country == 'CSK'] <- 'Eastern Europe'
border_tiles$region[border_tiles$country == 'DDR'] <- 'Eastern Europe'
border_tiles$region[border_tiles$country == 'YDYE'] <- 'Western Asia'
border_tiles$region[border_tiles$country == 'YUG'] <- 'Southern Europe'


border_tiles <- border_tiles %>% group_by(region, year) %>%
  dplyr::summarise(region_sentiment = sum(mean_predicted_sentiment, na.rm=TRUE), 
                   region_anxiety = sum(mean_anxiety, na.rm=TRUE),
                   region_anger = sum(mean_anger, na.rm=TRUE),
                   n_region = sum(!is.na(mean_predicted_sentiment))) %>%
  right_join(border_tiles, by=c('region', 'year')) %>% 
  mutate(region_sentiment = (region_sentiment - mean_predicted_sentiment)/n_region,
         region_anxiety = (region_anxiety - mean_anxiety)/n_region,
         region_anger = (region_anger - mean_anger)/n_region)

# world sentiment
border_tiles <- border_tiles %>% group_by(year) %>% dplyr::summarise(global_sentiment = mean(mean_predicted_sentiment, na.rm=TRUE)) %>%
  right_join(border_tiles, by = 'year')

# fig 10, sentiment, w/ year fixed effects and w/ neighbor sentiment
lfe_mod_sentiment <- felm(mean_predicted_sentiment ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                            Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                            I(10*Economic_Globalization) + I(10*Political_Globalization) +
                            I(Net_Migration_Prop_Imputed*100)*Polity + region_sentiment | country + year | 0 | country + year, data = border_tiles)

s <- summary(lfe_mod_sentiment)$coefficients
coeffs <- data.frame(s)
coeffs$varname <- rownames(s)
coeffs$lag <- '0-year lag'

border_tiles <- border_tiles[order(border_tiles$year),]
border_tiles <- border_tiles[order(border_tiles$country),]

# lagged versions
for(k in 1:3){
  border_tiles_lagged <- border_tiles %>% select(Interstate_Disputes,
                                                 Civil_Disputes_Neighbor,
                                                 Log_Population,
                                                 Log_GDP,
                                                 Economic_Globalization,
                                                 Political_Globalization,
                                                 Net_Migration_Prop_Imputed,
                                                 Polity,
                                                 region_sentiment, 
                                                 mean_predicted_sentiment,
                                                 year,
                                                 country) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, mean_predicted_sentiment), 
                         .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_sentiment <- felm(mean_predicted_sentiment ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                              Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                              I(10*Economic_Globalization) + I(10*Political_Globalization) +
                              I(Net_Migration_Prop_Imputed*100)*Polity + region_sentiment | country + year | 0 | country + year, 
                            data = border_tiles_lagged)
  
  s <- summary(lfe_mod_sentiment)$coefficients
  coeffs_lag <- data.frame(s)
  coeffs_lag$varname <- rownames(s)
  coeffs_lag$lag <- paste(k, '-year', sep='')
  
  coeffs <- rbind(coeffs,
                  coeffs_lag)
}

coeffs <- coeffs %>% mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, lag) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, 
                          `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          mean_predicted_sentiment = 'Sentiment',
                          `I(Net_Migration_Prop_Imputed * 100)` = 'Net Migration (%)',
                          `I(Net_Migration_Prop_Imputed * 100):Polity` = 'Net Migration (%):Polity',
                          region_sentiment = 'Regional Average Sentiment')) %>%
  mutate(varname = factor(varname, levels=rev(c('Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'Regional Average Sentiment',
                                                'log(Population)', 'log(GDP/Capita)')))) 

#Figure H2: Full results for in-text sentiment model (Text Figure 10)
gg <- ggplot(coeffs) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=lag), position = position_dodge(width=0.2), width=0.2) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.2), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Sentiment DV, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  theme(legend.title = element_blank())

ggsave(gg, filename = 'Figure H2.png', dpi=1000, width = 5, height=4, bg="white")

# fig 10, panel 1
gg_sentiment <- ggplot(coeffs %>% filter(!varname %in% c('log(Population)', 'log(GDP/Capita)'), lag=='1-year')) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper), position = position_dodge(width=0.1), width=0.1) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.1), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Sentiment DV, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  theme(legend.title = element_blank())

# fig 10 variant, anxiety, w/ year fixed effects and w/ neighbor anxiety
lfe_mod_anxiety <- felm(mean_anxiety ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                          Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                          I(10*Economic_Globalization) + I(10*Political_Globalization) +
                          I(Net_Migration_Prop_Imputed*100)*Polity + region_anxiety | country + year | 0 | country + year, data = border_tiles)

s <- summary(lfe_mod_anxiety)$coefficients
coeffs <- data.frame(s)
coeffs$varname <- rownames(s)
coeffs$lag <- '0-year lag'

border_tiles <- border_tiles[order(border_tiles$year),]
border_tiles <- border_tiles[order(border_tiles$country),]

# lagged versions
for(k in 1:3){
  border_tiles_lagged <- border_tiles %>% select(Interstate_Disputes,
                                                 Civil_Disputes_Neighbor,
                                                 Log_Population,
                                                 Log_GDP,
                                                 Economic_Globalization,
                                                 Political_Globalization,
                                                 Net_Migration_Prop_Imputed,
                                                 Polity,
                                                 region_anxiety, 
                                                 mean_anxiety,
                                                 year,
                                                 country) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, mean_anxiety), 
                         .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_anxiety <- felm(mean_anxiety ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                            Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                            I(10*Economic_Globalization) + I(10*Political_Globalization) +
                            I(Net_Migration_Prop_Imputed*100)*Polity + region_anxiety | country + year | 0 | country + year, 
                          data = border_tiles_lagged)
  
  
  s <- summary(lfe_mod_anxiety)$coefficients
  coeffs_lag <- data.frame(s)
  coeffs_lag$varname <- rownames(s)
  coeffs_lag$lag <- paste(k, '-year', sep='')
  
  coeffs <- rbind(coeffs,
                  coeffs_lag)
}

coeffs <- coeffs %>% mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, lag) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, 
                          `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          mean_anxiety = 'Anxiety',
                          `I(Net_Migration_Prop_Imputed * 100)` = 'Net Migration (%)',
                          `I(Net_Migration_Prop_Imputed * 100):Polity` = 'Net Migration (%):Polity',
                          region_anxiety = 'Regional Average Anxiety')) %>%
  mutate(varname = factor(varname, levels=rev(c('Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'Regional Average Anxiety',
                                                'log(Population)', 'log(GDP/Capita)')))) 

# Appendix figure J2 for mean_anxiety 
gg <- ggplot(coeffs) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=lag), position = position_dodge(width=0.2), width=0.2) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.2), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  guides(color = guide_legend(reverse = TRUE)) + 
  theme(legend.title = element_blank())

ggsave(gg, filename = 'Figure J2.png', dpi=1000, width = 5, height=4, bg="white")

# fig 10, panel 2 
gg_anxiety <- ggplot(coeffs %>% filter(!varname %in% c('log(Population)', 'log(GDP/Capita)'), lag=='1-year')) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper), position = position_dodge(width=0.1), width=0.1) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.1), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Sentiment DV, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  theme(legend.title = element_blank())

# fig 10, 2 panels
library(gridExtra)
grid.arrange(gg_sentiment, gg_anxiety, ncol=2)

# Appendix figure J1: sentiment, w/o year fixed effects and w/ world
lfe_mod_sentiment <- felm(mean_predicted_sentiment ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                            Log_Population + Log_GDP + #I(Population/10000000) + I(GDP/10000) + 
                            I(10*Economic_Globalization) + I(10*Political_Globalization) +
                            I(Net_Migration_Prop_Imputed*100)*Polity + world_sentiment | country | 0 | country, data = border_tiles)

s <- summary(lfe_mod_sentiment)$coefficients
coeffs <- data.frame(s)
coeffs$varname <- rownames(s)
coeffs$lag <- '0-year'

border_tiles <- border_tiles[order(border_tiles$year),]
border_tiles <- border_tiles[order(border_tiles$country),]

# lagged versions
for(k in 1:3){
  border_tiles_lagged <- border_tiles %>% select(Interstate_Disputes,
                                                 Civil_Disputes_Neighbor,
                                                 Log_Population,
                                                 Log_GDP,
                                                 Economic_Globalization,
                                                 Political_Globalization,
                                                 Net_Migration_Prop_Imputed,
                                                 Polity,
                                                 world_sentiment, 
                                                 mean_predicted_sentiment,
                                                 year,
                                                 country) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, mean_predicted_sentiment), 
                         .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_sentiment <- felm(mean_predicted_sentiment ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                              Log_Population + Log_GDP + 
                              I(10*Economic_Globalization) + I(10*Political_Globalization) +
                              I(Net_Migration_Prop_Imputed*100)*Polity + world_sentiment | country | 0 | country, 
                            data = border_tiles_lagged)
  
  
  s <- summary(lfe_mod_sentiment)$coefficients
  coeffs_lag <- data.frame(s)
  coeffs_lag$varname <- rownames(s)
  coeffs_lag$lag <- paste(k, '-year', sep='')
  
  coeffs <- rbind(coeffs,
                  coeffs_lag)
}

coeffs <- coeffs %>% mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, lag) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, 
                          `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          mean_predicted_sentiment = 'Sentiment',
                          `I(Net_Migration_Prop_Imputed * 100)` = 'Net Migration (%)',
                          `I(Net_Migration_Prop_Imputed * 100):Polity` = 'Net Migration (%):Polity',
                          world_sentiment = 'World Average Sentiment')) %>%
  mutate(varname = factor(varname, levels=rev(c('Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'World Average Sentiment',
                                                'log(Population)', 'log(GDP/Capita)')))) 

gg <- ggplot(coeffs) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=lag), position = position_dodge(width=0.2), width=0.2) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.2), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  guides(color = guide_legend(reverse = TRUE),
         name='Lags') + 
  theme(legend.title = element_blank())

ggsave(gg, filename = 'Figure J1.png', dpi=1000, width = 5, height=4, bg="white")

# Figure J3: fig 10 variant, anger, w/ year fixed effects and w/ neighbor anger
lfe_mod_anger <- felm(mean_anger ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                        I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                        I(Net_Migration_Prop_Imputed*100)*Polity + region_anger | country + year | 0 | country + year, data = border_tiles)
s <- summary(lfe_mod_anger)$coefficients
coeffs <- data.frame(s)
coeffs$varname <- rownames(s)
coeffs$lag <- '0-year lag'

border_tiles <- border_tiles[order(border_tiles$year),]
border_tiles <- border_tiles[order(border_tiles$country),]

# lagged versions
for(k in 1:3){
  border_tiles_lagged <- border_tiles %>% select(Interstate_Disputes,
                                                 Civil_Disputes_Neighbor,
                                                 Log_Population,
                                                 Log_GDP,
                                                 Economic_Globalization,
                                                 Political_Globalization,
                                                 Net_Migration_Prop_Imputed,
                                                 Polity,
                                                 region_anger, 
                                                 mean_anger,
                                                 year,
                                                 country) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, mean_anger), 
                         .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_anger <- felm(mean_anger ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                          I(10*Economic_Globalization) + I(10*Political_Globalization) +
                          I(Net_Migration_Prop_Imputed*100)*Polity + region_anger | country + year | 0 | country + year, 
                        data = border_tiles_lagged)
  
  
  s <- summary(lfe_mod_anger)$coefficients
  coeffs_lag <- data.frame(s)
  coeffs_lag$varname <- rownames(s)
  coeffs_lag$lag <- paste(k, '-year', sep='')
  
  coeffs <- rbind(coeffs,
                  coeffs_lag)
}

coeffs <- coeffs %>% mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, lag) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, 
                          `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          mean_anger = 'Anger',
                          `I(Net_Migration_Prop_Imputed * 100)` = 'Net Migration (%)',
                          `I(Net_Migration_Prop_Imputed * 100):Polity` = 'Net Migration (%):Polity',
                          region_anger = 'Regional Average Anger')) %>%
  mutate(varname = factor(varname, levels=rev(c('Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'Regional Average Anger'))))

ggplot(coeffs) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=lag), position = position_dodge(width=0.2), width=0.2) + 
  geom_point(aes(y=varname, x=Estimate, group = lag), position=position_dodge(width=0.2), size=1) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Sentiment DV, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  theme(legend.title = element_blank())

#################
# walls as a DV
#################
walls <- read_csv('Sentiments_border_walls_2021.05.07.csv')

cown_iso3 <- c('260' = 'DEU', '265' = 'DDR', '315' = 'CSK', '345' = 'YUG', 
               '347' = 'UNK', '678' = 'YEM', '680' = 'YDM', '816' = 'DVM')

walls <- walls %>% mutate(state_a = countrycode(as.character(ccode_a), 'cown', 'iso3c', nomatch=NULL),
                          state_b = countrycode(as.character(ccode_b), 'cown', 'iso3c', nomatch=NULL))

walls$state_a <- dplyr::recode(walls$state_a, !!!cown_iso3)
walls$state_b <- dplyr::recode(walls$state_b, !!!cown_iso3)

walls <- walls %>% filter(wall == 1)
walls$construction_end <- replace_na(walls$construction_end, 2017)
walls$construction_start <- replace_na(walls$construction_start, 1970)

walls <- walls %>% rowwise() %>%
  do(data.frame(country_1 = .$state_a, country_2 = .$state_b,
                year= seq(.$construction_start, .$construction_end), wall = .$wall))

walls <- walls %>% group_by(country_1, year) %>% dplyr::summarise(walls = sum(wall))

border_tiles <- left_join(border_tiles, walls, by=c('country' = 'country_1', 'year'))
border_tiles$walls <- replace_na(border_tiles$walls, 0)
border_tiles$walls <- ifelse(border_tiles$walls > 0, 1, 0)

border_tiles_joined <- border_tiles %>% group_by(year) %>% 
  dplyr::summarise(world_sentiment = mean(mean_predicted_sentiment, na.rm=TRUE)) %>%
  right_join(border_tiles, by='year')

library(mice)
# original version
coeffs_sentiment <- data.frame()
lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                        Log_GDP +  
                        I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                        I(100*Net_Migration_Prop_Imputed)*Polity + 
                        mean_predicted_sentiment | country + year | 0 | country + year, 
                      data = border_tiles_joined)

s <- data.frame(summary(lfe_mod_walls)$coefficients)
s$model <- '0-year'
s$varname <- rownames(s)
coeffs_sentiment <- rbind(coeffs_sentiment, s)

border_tiles_joined <- border_tiles_joined[order(border_tiles_joined$country, border_tiles_joined$year), ]

# lagged versions - mean_sentiment
for(k in 1:3){
  border_tiles_lagged <- border_tiles_joined %>% select(Interstate_Disputes,
                                                        Civil_Disputes_Neighbor,
                                                        Log_Population,
                                                        Log_GDP,
                                                        Economic_Globalization,
                                                        Political_Globalization,
                                                        Net_Migration_Prop_Imputed,
                                                        Polity,
                                                        mean_predicted_sentiment,
                                                        mention_binary_undergrad, 
                                                        year, country, walls, #a number of emotions deleted - SL
                                                        mean_anxiety) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, walls), .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                          Log_GDP +  
                          I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                          I(100*Net_Migration_Prop_Imputed)*Polity + 
                          mean_predicted_sentiment  | country + year | 0 | country + year, 
                        data = border_tiles_lagged)
  
  s <- data.frame(summary(lfe_mod_walls)$coefficients)
  s$model <- glue::glue('{k}-year')
  s$varname <- rownames(s)
  coeffs_sentiment <- rbind(coeffs_sentiment, s)
  
}

coeffs_sentiment <- coeffs_sentiment %>% 
  mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, model) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          mean_predicted_sentiment = 'Sentiment',
                          `I(100 * Net_Migration_Prop_Imputed)` = 'Net Migration (%)',
                          `I(100 * Net_Migration_Prop_Imputed):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(varname = factor(varname, levels=rev(c('Sentiment', 
                                                'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'log(Population)', 'log(GDP/Capita)'))))

gg_walls_1 <- ggplot(coeffs_sentiment) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=model), width=0.5, position = position_dodge(width=0.2)) + 
  geom_point(aes(y=varname, x=Estimate, group=model), position = position_dodge(width=0.2), size=0.5) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Wall construction, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  scale_color_discrete(name = 'Lags') +
  theme(legend.position = 'none')

#mean_anxiety
coeffs_anxiety <- data.frame()
lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                        Log_GDP +  
                        I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                        I(100*Net_Migration_Prop_Imputed)*Polity + 
                        I(10*mean_anxiety) | country + year | 0 | country + year, 
                      data = border_tiles_joined)

s <- data.frame(summary(lfe_mod_walls)$coefficients)
s$model <- '0-year'
s$varname <- rownames(s)
coeffs_anxiety <- rbind(coeffs_anxiety, s)

border_tiles_joined <- border_tiles_joined[order(border_tiles_joined$country, border_tiles_joined$year), ]

# lagged versions - mean_anxiety
for(k in 1:3){
  border_tiles_lagged <- border_tiles_joined %>% select(Interstate_Disputes,
                                                        Civil_Disputes_Neighbor,
                                                        Log_Population,
                                                        Log_GDP,
                                                        Economic_Globalization,
                                                        Political_Globalization,
                                                        Net_Migration_Prop_Imputed,
                                                        Polity,
                                                        mean_predicted_sentiment,
                                                        mention_binary_undergrad, 
                                                        year, country, walls,  
                                                        mean_anxiety) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, walls), .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                          Log_GDP +  
                          I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                          I(100*Net_Migration_Prop_Imputed)*Polity + 
                          I(10*mean_anxiety) | country + year | 0 | country + year, 
                        data = border_tiles_lagged)
  
  s <- data.frame(summary(lfe_mod_walls)$coefficients)
  s$model <- glue::glue('{k}-year')
  s$varname <- rownames(s)
  coeffs_anxiety <- rbind(coeffs_anxiety, s)
  
}

coeffs_anxiety <- coeffs_anxiety %>% 
  mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, model) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          `I(10 * mean_anxiety)` = 'Anxiety',                          
                          `I(100 * Net_Migration_Prop_Imputed)` = 'Net Migration (%)',
                          `I(100 * Net_Migration_Prop_Imputed):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(varname = factor(varname, levels=rev(c('Anxiety', 
                                                'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'log(Population)', 'log(GDP/Capita)'))))

gg_walls_2 <- ggplot(coeffs_anxiety) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=model), width=0.5, position = position_dodge(width=0.2)) + 
  geom_point(aes(y=varname, x=Estimate, group=model), position = position_dodge(width=0.2), size=0.5) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Wall construction, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  scale_color_discrete(name = 'Lags') + theme(legend.position = 'right')

# figure J4, appendix
grid.arrange(gg_walls_1, gg_walls_2, nrow=1)

# fig 11, main text
gg_walls_1_main <- ggplot(coeffs_sentiment %>% filter(model == '1-year')) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper), width=0.1) + 
  geom_point(aes(y=varname, x=Estimate, group=model), size=0.5) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Wall construction, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  scale_color_discrete(name = 'Lags') +
  theme(legend.position = 'none')

gg_walls_2_main <- ggplot(coeffs_anxiety %>% filter(model == '1-year')) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper), width=0.1) + 
  geom_point(aes(y=varname, x=Estimate, group=model), size=0.5) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Wall construction, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  scale_color_discrete(name = 'Lags') + theme(legend.position = 'right')

grid.arrange(gg_walls_1_main, gg_walls_2_main, nrow=1)

# Figure J5
# imputed, sentiment, 1-year lagged
border_tiles_joined <- border_tiles_joined[order(border_tiles_joined$country, border_tiles_joined$year), ]
border_tiles_lagged <- border_tiles_joined %>% select(Interstate_Disputes,
                                                      Civil_Disputes_Neighbor,
                                                      Log_Population,
                                                      Log_GDP,
                                                      Economic_Globalization,
                                                      Political_Globalization,
                                                      Net_Migration_Prop_Imputed,
                                                      Polity,
                                                      mean_predicted_sentiment,
                                                      mention_binary_undergrad, 
                                                      year, country, walls, mean_anxiety) %>%
  ungroup() %>%
  dplyr::mutate(across(!c(year, country, walls), .fns = ~panel_lag(., id = country, time = year, lag=1)))
impute <- mice(border_tiles_lagged %>% select(Interstate_Disputes,
                                              Civil_Disputes_Neighbor,
                                              Log_Population,
                                              Log_GDP,
                                              Economic_Globalization,
                                              Political_Globalization,
                                              Net_Migration_Prop_Imputed,
                                              Polity,
                                              mean_predicted_sentiment,
                                              mention_binary_undergrad, 
                                              year, country, walls, mean_anxiety))
lfe_mod_walls <- with(impute, felm(walls ~  Interstate_Disputes + Civil_Disputes_Neighbor + 
                                     Log_GDP +  
                                     I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                                     I(100*Net_Migration_Prop_Imputed)*Polity + 
                                     I(10*mean_predicted_sentiment)  | country + year | 0 | country + year))


s <- data.frame(summary(pool(lfe_mod_walls)))

s_sentiment <- s %>% mutate(lower = estimate - 1.96*std.error, upper = estimate + 1.96*std.error) %>% 
  select(term, estimate, lower, upper, std.error) %>% mutate(term = recode(term, `I(10 * Political_Globalization)` = 'Political Globalization',
                                                                           `I(10 * Economic_Globalization)` = 'Economic Globalization',
                                                                           Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                                                                           Interstate_Disputes = 'Interstate Disputes',
                                                                           Log_Population = 'log(Population)',
                                                                           Log_GDP = 'log(GDP/Capita)',
                                                                           `I(10 * mean_predicted_sentiment)` = 'Sentiment',
                                                                           `I(100 * Net_Migration_Prop_Imputed)` = 'Net Migration (%)',
                                                                           `I(100 * Net_Migration_Prop_Imputed):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(term = factor(term, levels=rev(c('Sentiment', 
                                          'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                          'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                          'Economic Globalization', 'Political Globalization', 
                                          'log(Population)', 'log(GDP/Capita)'))))

sentiment_imputed <- ggplot(s_sentiment) + 
  geom_errorbar(aes(x=term, ymin=lower, ymax=upper), width=0.2) + 
  geom_point(aes(x=term, y=estimate)) + 
  geom_hline(aes(yintercept=0), linetype='dashed') + 
  theme_minimal() + coord_flip() + xlab(NULL) + ylab(NULL) 

# imputed, anxiety, 1-year lagged data
impute <- mice(border_tiles_lagged %>% select(Interstate_Disputes,
                                              Civil_Disputes_Neighbor,
                                              Log_Population,
                                              Log_GDP,
                                              Economic_Globalization,
                                              Political_Globalization,
                                              Net_Migration_Prop_Imputed,
                                              Polity,
                                              mean_predicted_sentiment,
                                              mention_binary_undergrad, 
                                              year, country, walls, mean_anxiety))
lfe_mod_walls <- with(impute, felm(walls ~  Interstate_Disputes + Civil_Disputes_Neighbor + 
                                     Log_GDP +  
                                     I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                                     I(100*Net_Migration_Prop_Imputed)*Polity + 
                                     I(10*mean_anxiety)  | country + year | 0 | country + year))


s <- data.frame(summary(pool(lfe_mod_walls)))

s_anxiety <- s %>% mutate(lower = estimate - 1.96*std.error, upper = estimate + 1.96*std.error) %>% 
  select(term, estimate, lower, upper, std.error) %>% mutate(term = recode(term, `I(10 * Political_Globalization)` = 'Political Globalization',
                                                                           `I(10 * Economic_Globalization)` = 'Economic Globalization',
                                                                           Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                                                                           Interstate_Disputes = 'Interstate Disputes',
                                                                           Log_Population = 'log(Population)',
                                                                           Log_GDP = 'log(GDP/Capita)',
                                                                           `I(10 * mean_anxiety)` = 'Anxiety',
                                                                           `I(100 * Net_Migration_Prop_Imputed)` = 'Net Migration (%)',
                                                                           `I(100 * Net_Migration_Prop_Imputed):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(term = factor(term, levels=rev(c('Anxiety', 
                                          'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                          'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                          'Economic Globalization', 'Political Globalization', 
                                          'log(Population)', 'log(GDP/Capita)'))))

anxiety_imputed <- ggplot(s_anxiety) + 
  geom_errorbar(aes(x=term, ymin=lower, ymax=upper), width=0.2) + 
  geom_point(aes(x=term, y=estimate)) + 
  geom_hline(aes(yintercept=0), linetype='dashed') + 
  theme_minimal() + coord_flip() + xlab(NULL) + ylab(NULL) 

grid.arrange(sentiment_imputed, anxiety_imputed, nrow=1)


# Figure J6: adding anger (same as anxiety)
# adding anger (same as anxiety)
coeffs_anger <- data.frame()
lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                        Log_GDP +  
                        I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                        I(100*Net_Migration_Prop_Imputed)*Polity + 
                        I(10*mean_anger)  | country + year | 0 | country + year, 
                      data = border_tiles_joined)

s <- data.frame(summary(lfe_mod_walls)$coefficients)
s$model <- '0-year'
s$varname <- rownames(s)
coeffs_anger <- rbind(coeffs_anger, s)

border_tiles_joined <- border_tiles_joined[order(border_tiles_joined$year),]
border_tiles_joined <- border_tiles_joined[order(border_tiles_joined$country),]

# lagged versions
for(k in 1:3){
  border_tiles_lagged <- border_tiles_joined %>% select(Interstate_Disputes,
                                                        Civil_Disputes_Neighbor,
                                                        Log_Population,
                                                        Log_GDP,
                                                        Economic_Globalization,
                                                        Political_Globalization,
                                                        Net_Migration_Prop_Imputed,
                                                        Polity,
                                                        mean_predicted_sentiment,
                                                        mention_binary_undergrad, 
                                                        year, country, walls,
                                                        mean_anger, mean_anxiety) %>%
    ungroup() %>%
    dplyr::mutate(across(!c(year, country, walls), .fns = ~panel_lag(., id = country, time = year, lag=k)))
  
  lfe_mod_walls <- felm(walls ~ Interstate_Disputes + Civil_Disputes_Neighbor + 
                          Log_GDP +  
                          I(10*Economic_Globalization) + I(10*Political_Globalization) + 
                          I(100*Net_Migration_Prop_Imputed)*Polity + 
                          I(10*mean_anger)  | country + year | 0 | country + year, 
                        data = border_tiles_lagged)
  
  s <- data.frame(summary(lfe_mod_walls)$coefficients)
  s$model <- glue::glue('{k}-year')
  s$varname <- rownames(s)
  coeffs_anger <- rbind(coeffs_anger, s)
  
}

coeffs_anger <- coeffs_anger %>% 
  mutate(lower = Estimate - 1.96*`Cluster.s.e.`, upper = Estimate + 1.96*`Cluster.s.e.`) %>% 
  select(varname, Estimate, lower, upper, `Cluster.s.e.`, model) %>%
  mutate(Estimate = Estimate, lower = lower, upper = upper) %>%
  mutate(varname = recode(varname, `I(10 * Political_Globalization)` = 'Political Globalization',
                          `I(10 * Economic_Globalization)` = 'Economic Globalization',
                          Civil_Disputes_Neighbor = 'Civil Disputes (Neighbor)',
                          Interstate_Disputes = 'Interstate Disputes',
                          Log_Population = 'log(Population)',
                          Log_GDP = 'log(GDP/Capita)',
                          `I(10 * mean_anger)` = 'Anger',
                          `I(100 * Net_Migration_Prop_Imputed)` = 'Net Migration (%)',
                          `I(100 * Net_Migration_Prop_Imputed):Polity` = 'Net Migration (%):Polity')) %>%
  mutate(varname = factor(varname, levels=rev(c('Anger', 
                                                'Interstate Disputes', 'Civil Disputes (Neighbor)', 
                                                'Polity', 'Net Migration (%)', 'Net Migration (%):Polity',
                                                'Economic Globalization', 'Political Globalization', 
                                                'log(Population)', 'log(GDP/Capita)'))))

gg_walls_3 <- ggplot(coeffs_anger) + 
  geom_errorbar(aes(y=varname, xmin=lower, xmax=upper, color=model), width=0.5, position = position_dodge(width=0.2)) + 
  geom_point(aes(y=varname, x=Estimate, group=model), position = position_dodge(width=0.2), size=0.5) + 
  geom_vline(aes(xintercept=0), linetype='dashed') + 
  theme_minimal() + xlab(NULL) + ylab(NULL) + 
  #ggtitle('Wall construction, country-year (linear probability model, SE scale)') + 
  guides(color = guide_legend(reverse = TRUE)) + 
  scale_color_discrete(name = 'Lags') + theme(legend.position = 'right')

grid.arrange(gg_walls_1, gg_walls_3, nrow=1)

